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We show the results of two dimensional general relativistic inviscid and isothermal hydrodynamical simulations comparing the 
behavior of co-rotating (with respect to the black hole rotation) and counter-rotating circumbinary quasi-Keplerian discs in the post 
merger phase of a supermassive binary black hole system. While confirming the spiral shock generation within the disc due to the 
combined effects of mass loss and recoil velocity of the black hole, we find that the maximum luminosity of counter-rotating discs 
is a factor ~ (2 - 12) higher than in the co-rotating case, depending on the spin of the black hole. On the other hand, the luminosity 
peak happens ~ 10 days later with respect to the co-rotating case. Although the global dynamics of counter-rotating discs in the 
post merger phase of a merging event is very similar to that for co-rotating discs, an important difference has been found. In fact, 
increasing the spin of the central black hole produces more luminous co-rotating discs while less luminous counter-rotating ones. 
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1. Introduction 
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> In the last few years a lot of attention has been given to the 
^ possibility of detecting the electromagnetic (EM) counterpart of 
I— ' the gravitational wave signal produced in the merger event be- 
tween two supermassive binary black holes (SMBBHs), whose 
gravitational waves should be detected by the planned Laser 
'Interferometric Space Antenna (LISA) (see, among the others . 
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Bogdano vic et al.L l201 1: Haiman et al., 2009; Palenz uela et al 
20ia 



i .Sesana et al.l {1009l .Stavridis et al.. ,2009; .Babak et al 



120111) . 



Our current understanding of the physics of these systems 
distinguishes between a pre-merger phase and a post-merger 
phase. The EM emission properties of the pre- merger phase 
have been studied by several authors, includ i ng Mosta et al.J 
bOlOhi'Bode et aL'^OlOh Jchang et alJ(l201(il):IShaDi ro (2010): 

in particular. 



(2010). IChang et al 



(l2010l) 



showed that, approximately one day before the merger is com- 
pleted, there is a late time precursor brightening from tidal and 
viscous dissipation in the inner disc. On the other hand, the typ- 
ical behavior of the SMBBH system in the post-merger phase 
is characterized by two distinct but related physical effects hav- 
ing to do with the emission of gravitational waves, namely an 
abrupt mass loss and a sudden recoil of the resultant black hole. 
Observationally, the "recoil event" has recently found indirect 
confirmation throug h spectropolarirnetric m easurements of the 



quasar E1821+643 jRobinson et al.ll2010l) . that are consistent 



with a scattering model in which the broad-line region has two 
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components, moving with different bulk velocities away from 
the observer. The kick velocity and the mass loss of the black 
hole are extremely relevant both because they represent an im- 
print of the gravitational wave emission, and because they have 
strong dynamical effects on the circumbinary disc, whose inner 
radius is that possessed at the time when the binary in-spiral 
time scale became shorter than the viscous time scale in the 
disc. 

When investigating the response of the disc to the mass loss 
of the black hole and its recoil velocity, an additional physi- 
cal effect that needs to be taken into account is that produced 
by counter-rotation. The term co-rotation (counter-rotation) 
has two different meanings when applied to the pre-merger 
phase, when it has to be intended as rotation in the same 
(opposite) sense of the merging binary, and when applied to 
the post-merger phase, when it has to be intended as rotation 
in the same (opposite) sense of the spinning resultant black 
hole. It is worth recalling that, on a broader sense, counter- 
rotat ing discs are a relativ ely common phenomenon in galax- 
ies. ICoccato et al. ( 2011 ) report a list of galaxies, in cluding 
the well studied case of NGC4550 ('R ubin et al.lll992h . where 
two co-spatial stellar discs, one orbiting prograde, the other 
one orbiting retrograde, have been observed. In general, large 
scale counter-rotating stellar discs manifest preferentially in 
early-type spirals, while 20% of the gas discs in SO galax- 
ies is thought to counter-rotate. As far as SMBBH s systems 
are concerned. iMacFadven and M ilosavlievic (2008) speculate 
about the important dynamical differences that may be present 
in a cou nter-rotating disc d uring the pre-merger phase. Quite 
recentlv [Nixon et al (12010^ have performed a detailed analy- 
sis, both analytical and numerical, of the retrograde accretion 
in merging supermassive black holes, reaching the conclusion 
that the binary orbit, by absorbing negative angular momentum 
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from the circumbinary disc via gas capture, may become ec- 
centric. In particular, the binary coalesces once it has absorbed 
the angular momentum of a gas mass comparable to that of the 
secondary black hole. 

All of these investigations, however, refer to the pre-merger 
phase. During the post-merger phase, it is not clear whether 
a retrograde disc may have significant effects on both the dy- 
namics and the EM signal emitted by the system. Several 
time dependent numerical investigations of prograde circumbi- 
nary discs in the po st-merger phase h a ve been performed . 



both in a New t onian jLippai et all 120081: lO'NeiU et all 120091: 



Corrales et all 2010; Rossi et al., 201C) and in a 



rel ativis- 

tic framework (|Megevand et al... .20091 : lAnderson et al. , 2010l : 
Zanotti et al. . 201ol) . With different emphasis and caveats, the 



majority of these works have confirmed the possibility of de- 
tecting a significant E M counterpart from SMBBHs. By ex- 



tending the analysis of IZanotti et al. I (l2010l) . in this paper the 



role of retrograde rotation in the post-merger phase is studied 
through numerical simulations, with the special aim of compar- 
ing the emitted luminosity with respect to the prograde rotation 
case. 

The paper is organized as follows. In Section|2]and in Sec. [3] 
the essential information about the physical properties of the 
initial models, and about the numerical code adopted in the sim- 
ulations are provided. Section|4]is devoted to the presentations 
of the results, while Section|5]contains a summary of the work. 
A signature {-, +, +, +} for the space-time metric is assumed 
and Greek letters (running from to 3) for four-dimensional 
space-time tensor components are used, while Latin letters 
(running from 1 to 3) will be employed for three-dimensional 
spatial tensor components. Moreover, c = G - 1 and the geo- 
metric system of units is extended by setting nip I kg - 1, where 
nip is the mass of the proton, while kg is the Boltzmann con- 
stant. In this way the temperature is a dimensionless quantity. 



counter-rotating, while maintaining all of the other physical pa- 
rameters unmodified. In particular, a specific counter-rotating 
model is larger than the corresponding co-rotating model with 
the same disc temperature. Another relevant difference among 
the two classes of models, which is clearly visible by inspec- 
tion of Table [T] is that, while for co-rotating models the radius 
of the maximum rest mass density decreases for larger spins of 
the central black hole, the opposite happens for counter-rotating 
model^. This will have important implications on the depen- 
dence of the emitted luminosity on the spin, as discussed in 
Sec.lH 

The recoil effect on the black hole resulting from a 



SMBBHs merger was predi cted well before jBekensteinl 11973 



Redmount and Rees , Il989h numerical relativistic simulations 
were able to prove its existence and measure its magnitude. Be- 
cause in this contribution only two dimensional simulations are 
performecd, the recoil velocity vector has to lie on the orbital 
plane of the disc, which is consistent with a binary merger be- 
tween two black holes having the same mass and with spins 
that are equal in mo dulus and both perpendicular to the orbita l 
angular momentum (ICampanelli et al.L 120071 : iRezzoUaL l2009l) . 
Since the recoil velocities in th e orbital plane are expected to 
be < 450km/s dKoppitz et all 120071: iHerrmann etall 12007 : 
Pollnev et al. , 20071 ). the conservative value Vk = 300km/s is 
adopted. In order to simulate this effect numerically, at time 
f = a Lorentz boost to the fluid velocity of the disc is applied, 
oriented along the radial direction with = 0. In addition, 2% 
of the total mass-energy of the merged black hole is assumed to 
be radiated away in gravitational waves. In practice, the initial 
model is first computed in the gravitational potential of the full 
black hole mass, and then evolved in the gravitational potential 
of the reduced mass. 



3. Numerical method 



2. Initial conditions 

The initial model is given by a stationary and axisymmetric 
configuration obtained after solving the relativistic Euler equa- 
tions in the fixed bac kground spacetirn e of a Kerr black hole. 
Such solution is due to Kozlowski et al. 1(11978') while a detailed 



description for co mputing the correspondin g equilibrium mod- 
els can be found in Daigne and FontI ( 2004 ). The specific angu- 
lar momentum of the resulting geometrically thick disc obeys a 
power law on the equatorial plane, namely i{r, - n/l) = Sr^, 
where S is positive or negative, corresponding to a disc rota- 
tion that is prograde or retrograde with respect to the black hole 
rotation. Table [T] reports the main parameters of the models 
considered, differing for the rotation law, for the spin of the 
central black hole, and for the temperature of the disc. Having 
chosen the exponent q to be close to 1/2 in all of the models, 
the rotation law tends to the Keplerian one, and the disc flattens 
towards the equatorial plane, thus making its vertical structure 
negligible. The equation of state of the initial model is that of a 
polytrope p = /^p^, with y = 5/3. It should be noted that it is not 
possible to build two models, one co-rotating and the other one 



The equations of general relativistic hydrodynamics are 
solved in the stationary spacetime of a Kerr black hole, 
writt en in Boyer-Lindquist coordinates, through the ECHO 
code jPel Zanna et alll2007l) . ECHO adopts a 3 + 1 split of space- 
time in which the space-time metric is decomposed according 
to 

ds^ = gpydx^dx" = -aW + jij (dx' + /3'dr)(dx-' + /?Mf), (1) 

where a is the lapse function,/?' is the shift vector, and y,j is the 
spatial metric tensor, with / and j running from 1 to 3. More- 
over, ECHO adopts a conservative formulation of the general- 
relativistic inviscid hydrodynamical equations, by solving the 
system 



(2) 



'Note that, when changing the spin of the blaclc hole, the potential gap be- 
tween the position of the cusp and the position of the inner radius is kept con- 



stant 

2 



It is worth recalling that the expansion of the disc along the vertical direc- 
tion is a very small effect, since vertical tangential velocities along the shock 
front are at least one order of magnitude smaller than those across the shock. 



2 



Table 1: Main properties of tlie initial models. From left to right the columns report the name of the model, the black hole spin parameter a, the power-law index q 
and the parameter .S of the angular momentum distribution, the inner and the outer radius of the disc, ri„ and Taut, the radius of the maximum rest-mass density r^, 
the orbital period at the radius of maximum rest-mass density t^, the uniform temperature T assumed during isothermal evolution. The adiabatic index y = 5/3 and 
the mass of the black hole is M = 1.0 X IO^Mq in each model. 
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where the vector of conservative variables 'ZY and the corre- 
sponding fluxes "F' in the ; direction are respectively given by 
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whereas the sources, in any stationary background metric, can 
be written as 





^aW^'djjik + Sidj0 - Udja 



(4) 



The determinant of the spatial metric -\fy and of the global met- 
ric V~? related by sfy = yf^ja. The conservative vari- 
ables {D,S j, U) that are evolved by the numerical scheme are 
related to the rest-mass density p, to the thermal pressure p and 
to the fluid velocity v, by 



D=pW, 
Si=phW\, 
U = phW^ - p. 



(5) 
(6) 
(7) 



where W = (1 - v^)"'^^ is the Lorentz factor of the fluid with 
respect to the Eulerian observer associated to the 3-1-1 splitting 
of the spacetime, and 



Wjj = phW^ViVj + pjij 



(8) 



is the fully spatial projection of the energy-momentum tensor 
of the fluid. 



The formulation above is appropriate for numerical inte- 
gration via standard high-resolution shock- captur i ng (H RSC) 
methods developed for the Euler equations (iToroill999h . All 
the simulations cover a physical duration much shorter than the 
viscous time-scale, making the inviscid approximation a very 
good one. In the following, however, we consider an isother- 
mal evolution of the hydrodynamical equations. Therefore, the 
temperature of the disc is assumed to be constant in space and 
in time, and set to the value T reported in Table [1] As a result, 
there is no need to evolve the energy equation for U, since the 
energy can be computed directly from the temperature and the 
latter is constant by construction. However, an equation for the 
time evolution of the internal energy is actually solved 



-pdti VyW') - p^/Ca Vy«') ' 



(9) 



where e is the specific internal energy and p - pe{y - 1). The 
solution of Eq. (|9]) allows to obtain a good estimate of the emit- 
ted luminosity, assuming that all the changes in the tempera- 
ture produced by local compressions will be dissipated away 
as radiat ion. This idea, prop osed in the Newtonian frame- 
work by Corrales et alj ( 20101) . has been fi rst implemented in 
a general relativistic context by IZanotti et al . (2010). The lu- 
minosity is then computed in the following way: the time- 
update of Eq. (|9]) provides the value of the quantity pe at time 
jM+i _ ^ where Af is the time-step of the simulation. After 
that, a volume integration of the difference p[e - e(r)] is per- 
formed, where e{T) - hgT l[mp(y - 1)] is the value prescribed 
by the assumption of isothermal evolution. Finally, division by 
Af provides the luminosity. The specific internal energy is reset 
to e = €{T) after each time-step, so as to guarantee that the evo- 
lution is isothermal. In practice, Eq. (|9]) is evolved in time with 
the only aim of computing the diff'erence p[e - e(r)], which is 
then assumed to be radiated instantaneously. 
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Figure 1 : Rest mass density distribution in tlie co-rotating model CoRot_a .0.9 (left panel) and in the counter-rotating model CounterRot_a .0.9 (right panel) for 
a recoil velocity Vi^ = 300km/s at the time (after merger) when the luminosity reaches a maximum, cgs units are reported. 



The radial numerical grid is discretised by choosing A^, = 
1200 points, non-uniformly distributed from rmin - lOOM to 
'"max = 6000M, while outflow boundary conditions are adopted 
both at r^in and r,^^^. The azimuthal grid extends from to 2jt, 
with symmetric (i.e. periodic) boundary conditions adopted at 
<p = 0, while the number of angular grid points is A^^ = 200. All 
runs are performed with a Courant-Friedrichs-Lewy coefficient 
CFL - 0.5. We have verified through a series of tests that the 
choice (A^;., A^^) = 1200 x 200 suflices to find a converged solu- 
tion, and, in consequence, all models have been computed with 
this canonical resolution. 



4. Results 



As shown numerically for the first time by iLippai et al 
(120081), and then confirmed repeatedly by several au- 



thors 
20101: 



2009; Anderson etal, 2010; 
the combined effects of mass loss 



blNeill et al.', '2009'; 'Corrales et al.', '2010"; 'Rossi et al.', 
iMeg evand et al., 
Zanotti et al.[ I2OIOI) . 
from the central black hole and of the recoil velocity generate 
a spiral shock pattern that transports angular momentum 
outwards, and become responsible of an enhanced luminositjU. 

Figure [T] shows the iso-density curves of the rest mass den- 
sity in the two models CoRot_a .0.9 and Count erRot_a .0.9, 
plotted when the emitted luminosity reaches its maximum. Be- 
cause the recoil velocity is not large, neither of the two discs 
has penetrated into the central cavity at this time. However, 
the spiral structure is still clearly visible, producing regions 
of high compression. The role of this spiral pattern has been 
somewhat overestimated in previous analysis. Indeed, as shown 
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Figure 2: Luminosity computed through the isothermal evolution 
for the co-rotating model CoRot_a.0.9 and the counter-rotating one 
CounterRot_a. . 9. 



^^Note, how ever, that the spiral shock can form even in the absence of a 
recoil velocity to'Neill et a"ill2009h . 



by IZanotti et al. I (l2010h through the use a sophisticated shock 
detector, such spiral pattern is not always a true physical shock. 
In addition, the same pattern can be obtained even in the ab- 
sence of a mass loss or of a recoil velocity, being produced 
as the result of non-axisymmetric hydrodynamical instabilities 
taking place in the disc. Therefore, the spiral pattern cannot 
be regarded as an unambiguous imprint of a recoiling black 
hole. It is also worth stressing that, although the inner edge of 
the accretion disc is placed in a region where general relativis- 
tic effects may be regarded as negligible, during the subsequent 
evolution (later than the time shown in Fig.[TJ, the central cavity 
is filled with gas, thus making the relativistic calculation perti- 
nent. 

Figure |2] on the other hand, shows the electromagnetic lu- 
minosity computed according to the procedure described in 
Sec. [3] The negative contributions that are produced in re- 
gions experiencing rarefactions have been neglected. Because 
when accounted for these negative contributions typically yield 
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Figure 3: Dependence of tlie luminosity on the spin of tlie black hole for corotating (left panel) and for counter-rotating (right panel) models. The insets show the 
peak luminosity as a function of the black hole spin. 



values that are a factor ~ 2 - 5 smaller, the values in Fig. ^ 
Fig. ID s hould be taken as upp er limits to the emitted luminos- 
ity (see ICorrales et al. , 20 id). More details regar ding this ap- 
proximation can be found in Zanotti et al. ( 2010l) . where such 
an approach is compared to alternative ones for extracting a 
luminosity in the absence of a fully radiation hydrodynamics 
treatment. The two curves refer to the models CoRot_a.0.9 
and CounterRot^. . 9, having the same temperature T = 
2.68 X lO^K and the same polytropic constant, subject to a 
recoil velocity - 300km/s and with a mass loss after the 
merging event that amounts to 2% of the initial mass of the 
black hole. As i nitially put into evidence by the analysis of 
Daigne and FontI (2004), counter-rotating models have larger 
extension with respect to co-rotating ones and, in addition, have 
a larger radius of the maximum rest mass density point. These 
two facts are ultimately responsible for the different behavior 
of the light curves reported in Fig. |2l Because is larger for 
CouiiterRot_a. . 9 than for CoRot_a.0.9, the spiral shock 
that forms after the merger event takes longer to reach the re- 
gion of maximum rest-mass density, which is also the region of 
highest compressions where most of the luminosity is emitted. 
As a result, the time of the maximum luminosity is fpeak ~ 32 
days for model CoRot_a.0.9 and fpeak ~ 41 days for model 
CounterRot_a. . 9. Moreover, the maximum luminosity at 
the peak is L ~ 1.3 x 10'*'*erg/s for model CoRot_a. . 9 and 
L ~ 3.4 X lO'^erg/s for model CounterRot^. .9. In a sec- 
ond and tightly related series of simulations the dependence 
of the light curves on the spin of the central black hole has 
been investigated, while keeping the same temperature within 
the discs. When the spin of the black hole is increased, the ra- 
dius of the maximum rest-mass density moves towards the 
center for co-rotating models, while it moves towards larger 
distances for counter-rotating models (see Tab. [TJ. Moreover, 



Table 2: First column: black hole spin. Second column: ratio of the peak lumi- 
nosity between counter-rotating models and co-rotating models. Third column: 
difference (in days) between the times of the peak in counter-rotating and in 
co-rotating discs. 
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(^Counter ^Coi"ot)peak 
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0.1 
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11.6 


0.5 


5.12 


11.5 


0.9 


2.45 


9.5 


0.99 


2.13 


9.2 



when the spin of the black hole is increased, the size of the disc 
increases for co-rotating models, while it decreases for counter- 
rotating models. Because of this, changing the spin of the black 
hole has opposite effects on the emitted luminosity for the two 
classes of models, as reported in Fig. [3] Indeed, the fluid com- 
pression, and therefore the energy released, is larger if it takes 
place deeper in the potential well, where the effects of the spi- 
ral shock are stronger This explains why increasing the spin 
of the central black hole produces more luminous co-rotating 
discs and less luminous counter-rotating ones, as schematically 
represented by the insets of Fig. [3] Tab. |2] provides an addi- 
tional information, by showing the ratio of the peak luminosity 
between counter-rotating and co-rotating discs, in terms of the 
spin of the black hole. Because of the anti-correlation high- 
lighted above, the discrepancy in the luminosity between the 
two classes of models is higher at smaller spins. 

Finally, in a last series of simulations the dependence of the 
emitted luminosity on the disc temperature has been consid- 
ered. The results of these analysis have been reported in Fig.|4l 
which shows the light curves for five models having the same 



5 



45 



- 

- / /■ 


\ 


1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 

' ' T = 4 47x1 

T=8.94xl05 _ 

1=2.68x106 _ 

T=8. 05x106 

T=l. 61x10^ _ 










-// / 
4'/ 

' 1 ' ^— 








-/''/ 
-''/ 

; / 
/ / 








W 45 

- / J 

/ M 
O 








44 









1 1 1 1 1 




6 6.5 7 
logT 

1 1 1 1 1 1 1 1 1 1 1 1 





50 100 150 200 250 300 



time after merger (days) 

Figure 4: Light curves for different values of the disc temperature in counter- 
rotating models. 

rotation law and the same spin of the central black hole, while 
different disc temperatures. An almost linear increase of the 
peak luminosity (see inset of Fig.lUi with the disc temperature 
has been found, with Lpeak ~ 1.3x IQ-^^erg/s for T ~ 1.6x10^/:. 
Such a linear increase is not a surprising result, given the fact 
that the luminosity is computed from local variations of the spe- 
cific internal energy and that the latter scales like e oc r. Finally, 
the linear increase of the peak luminosity with the disc temper- 
ature, which has been reported in Fig. |4] for counter-rotating 
models, is not affected by the rotation law, and the same linear 
dependence has been found also for co-rotating models. 

5. Conclusions 

We have performed two dimensional (on the orbital plane) 
relativistic isothermal hydrodynamical simulations of counter- 
rotating circumbinary discs that react to the mass loss and to the 
recoil velocity of the central black hole after a merging event 
of supermassive black holes. These calculations are very rel- 
evant for correlating the electromagnetic signal to the gravita- 
tional one, in view of the planned Laser Interferometric Space 
Antenna (LISA) mission. In this respect, after considering 
the dynamical response of a disc to a typical recoil velocity 
Vk = 300km/s oriented on the orbital plane, the results show 
that the maximum luminosity of counter-rotating discs is a fac- 
tor ~ 2 higher than in co-rotating discs for high spinning black 
holes, while it is a factor ~ 12 higher for slowly rotating black 
holes. On the other hand, the peak of the luminosity typically 
happens ~ 10 days later with respect to the co-rotating case, 
with only a weak dependence on the spin of the black hole. All 
of these effects are due to the different sizes and characteristic 
radii in the two classes of models. In particular, when the spin 
of the resultant black hole is larger, co-rotating discs in the post- 



merger phase are more luminous, while counter-rotating discs 
are less luminous. 
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